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Abstract 


The spiders in the genus Pardosa C.L. Koch, 1847, are a young lineage of the family Lycosidae Sundevall, 1833, that exhibit high 
species diversity and widespread distribution. Pardosa is abundant in open and disturbed environments. In fact, most of its species 
live in grasslands, and the few that live in forests switched habitats relatively recently. The genus markedly prefers grasslands with 
a broad range of climates. Thus, its origin and diversification were probably associated with grassland expansions during the late 
Cenozoic. To test this hypothesis, we developed a global phylogenetic hypothesis that helps reconstruct the biogeographic patterns 
of the genus Pardosa using three nuclear (18S, ITS2, and H3) and four mitochondrial (12S, 16S, NADH1, and COJ) loci. Our phylo- 
genetic analyses cover 133 (125 described and 8 as yet undescribed) grassland species of Pardosa using Trochosa ruricola (De Geer, 
1778) and Lycosa coelestis L. Koch, 1878, as outgroups. The results show that our selection of species in the genus is divided into 
four major clades: Clade I includes only P. crassipalpis Purcell, 1903, from South Africa; Clade II consists of a north-east African 
group (2 species) and a south-east Asian group (21 species); Clade III contains only P. sutherlandi (Gravely, 1924) from SE Asia; 
and Clade IV includes five species groups from Asia, Europe, and the Americas. The spiders of the genus probably originated in 
southern Africa or southern and eastern (SE) Asia at the Middle Miocene Climatic Optimum, about 19.40—14.18 Ma, and then ex- 
panded northwards to North America via the Bering Strait, as well as southwards to north-east Africa via the Arabian Peninsula, and 
westwards to Europe via western Asia between about 10.59 and 5.28 Ma. At least three exchanges occurred between North America 
and SE Asia, and at least two between Europe and North America. The biogeography of Pardosa in the past 14.18 Ma, associated 
with the evolution of grasses, suggested a late Cenozoic diversification of the genus as grasslands expanded. 
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Introduction 


Plants in the family Poaceae, usually called grasses, are 
composed of more than 12 thousand species (Govaerts 
et al. 2021). They originated during the late Cretaceous 
(about 100 million years ago, Ma; Gallaher et al. 2022) 
and currently cover over 25% of all land on Earth, com- 
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prising 35% of the subtropics and tropics. Grasses also 
occur in all terrestrial habitats with a broad range of cli- 
mates, from cold to hot and arid to wet. Climate, soils, 
fire, and herbivory shaped global grassy systems, includ- 
ing the C3 and C4 photosynthesis groups (Linder et al. 
2018; Schubert et al. 2019; Stromberg and Staver 2022). 
The C3 grasses developed the tolerance needed to survive 
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in colder climates (Schubert et al. 2019); the C4 grasses, 
evolutionarily derived from the C3 grasses, can prosper 
in hot and dry areas (Linder et al. 2018). Grassy biomes 
(>20 Ma; Stromberg and Staver 2022), associated with 
open-canopy habitats, harbor many spiders of the family 
Lycosidae Sundevall, 1833 (wolf spiders; Piacentini and 
Ramirez 2019). They interact ecologically with a diverse 
part of the lycosid fauna, both above and below ground; 
genera with representatives occurring in grasslands in- 
clude Pardosa C. L. Koch, 1847; Lycosa Latreille, 1804, 
Trochosa C. L. Koch, 1847; and Hippasa Simon, 1885. 
These lycosids rely on native grasses as a substrate to live 
on and to forage on grass-eating insects. 

Based on the remarkable abundance of wolf spiders 
in open habitats and in the fossil record (Wunderlich 
2004a), the retention of the third tarsal claw, and the con- 
servative morphology, Jocqué and Alderweireldt (2005) 
proposed that Lycosidae had a relatively recent origin 
and had co-evolved with grasslands and, by extension, 
all types of open habitats. Recently, a dated phylogenetic 
analysis (Piacentini and Ramirez 2019) inferred that the 
family appeared before the expansion of grasslands and 
diversified (about 50 Ma) with the reduction of tropical 
forests and the advance of open habitats. However, the 
linkages between the temporal-spatial diversification of 
Lycosidae, especially its most diverse genus Pardosa 
(based on the number of described species), and grass- 
land expansions have not been specifically tested. 

Pardosa wolf spiders inhabit nearly all terrestrial 
habitats worldwide and currently comprise 532 species 
(World Spider Catalog 2024). They are free-roaming 
predators that can disperse over short distances; the spid- 
erlings of some species have been reported to disperse by 
ballooning on silken threads (e.g., Richter 1970; Green- 
stone 1982). They are abundant in grasslands and all 
types of open habitats, such as wetlands, stream banks, 
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floodplains, glades, open deserts, farmlands, and human 
settlements (Jocqué and Alderweireldt 2005). The marked 
preference for grasslands contributed to their diversifica- 
tion and successful expansions when the arid and grass- 
lands extended over the globe. Previous studies indicated 
that Pardosa diversified between about 14 and 10 Ma 
(Piacentini and Ramirez 2019), and most of its species 
are remarkably distributed in the four areas with different 
grassland environments and grassy (C3 and C4) evolu- 
tionary histories: (1) North America, NA; (2) eastern and 
southern Africa, EA; (3) central, eastern, and southern 
Asia, CA; and (4) the region surrounding the Mediterra- 
nean Sea, RM. (Fig. 1; Suppl. material 1: table S1; e.g., 
Vogel 1964, 1970, 2004; Hanggi et al. 1995; Ivanov et 
al. 2023). Therefore, Pardosa is a particularly appropriate 
animal group to test the hypothesis of a diversification of 
specific lycosid lineages as grasslands expanded. Howev- 
er, the genus lacks a comprehensive phylogenetic hypoth- 
esis, and thus its lineage compositions, phylogenetic re- 
lationships, and geographical distribution remain unclear. 

Spiders of the genus Pardosa have been proposed as 
effective biological control agents for pests in agricultur- 
al systems globally. As the biodiversity and ecological 
prominence of Pardosa spiders are increasingly recog- 
nized, there is a demand for illumination of their biogeo- 
graphical patterns and diversification mechanisms for the 
purposes of protection and ecological management. The 
purpose of this study was to produce a robust phyloge- 
netic hypothesis to elucidate the relationship of Pardosa 
grassland lineages based on a global sample of species 
using three nuclear (18S, ITS2, and H3) and four mito- 
chondrial (12S, 16S, NADHI/, and CO/) loci. Further- 
more, we explored the origin, diversification timeline, 
and global expansion history using a dated phylogenetic 
tree and tested the co-evolution between Pardosa spiders 
and grasslands. 


Figure 1. The major distributional areas of Pardosa spiders and the main biomes across the globe. White spots indicate the probable 
type localities of all known Pardosa species; red spots indicate the sample localities of Pardosa spiders sequenced in this study. 


zse.pensoft.net 


Zoosyst. Evol. 100 (4) 2024, 1287-1296 


Materials and methods 
Sampling and sequence analyses 


We collected 27 Pardosa grassland species across Chi- 
na between 2018 and 2024 (Fig. 1; Suppl. material 1: 
table S1). They were fixed in absolute ethanol and stored 
at -20 °C prior to analysis. Voucher specimens were de- 
posited at Shaoxing University, Shaoxing, China. TI- 
ANamp genomic DNA kits were used to extract DNA. 
Sequences of nuclear (18S, ITS2, and H3) and mitochon- 
drial (12S, 16S, COL, and NADH1) loci were amplified 
using seven pairs of primers (Suppl. material 1: table 
S2). Sequencing was performed utilizing BigDye tech- 
nology on an ABI 3730 automated sequencer (Applied 
Biosystems, Forter City, CA, USA). BioEdit (Hall 1999) 
and SEQUENCHER v4.1.2. (Gene Codes, USA) were 
applied to sequence quality checking and editing. New 
sequences were submitted to GenBank (Suppl. material 
1: table S1). Additional sequences of the Pardosa grass- 
land species from North America, Europe, East Asia, 
and Africa were taken from Zehethofer and Sturmbau- 
er (1998), Muster and Berendonk (2006), Murphy et al. 
(2006), Piacentini and Ramirez (2019), and Roslin et al. 
(2022), as well as GenBank and BOLDSYSTEM (Suppl. 
material 1: table S1). Alignments of genes were generat- 
ed using Clustal X (Jeanmougin et al. 1998). Sequences 
of the seven genes were merged to perform phylogenetic 
analyses, divergence time estimates, and biogeographical 
reconstructions. The data matrix will be submitted to the 
Dryad database (online at http://datadryad.org/). 


Phylogenetic analyses 


Phylogenetic relationships among the sampled Pardosa 
Species were inferred using both the maximum likelihood 
and Bayesian inference approaches. The wolf spiders 
Trochosa ruricola (De Geer, 1778) and Lycosa coelestis 
L. Koch, 1878, were used as outgroups. Maximum like- 
lihood analyses were implemented using the fast online 
phylogenetic tool W-IQ-TREE (Trifinopoulos et al. 2016). 
The optimal substitution model for each gene partition 
(TIM2+F+I+G4 for 12S, 16S, NADH], and COV, K2P+I 
for 18S; K2P+R2 for H3; and K2P+I for ITS2) was esti- 
mated simultaneously using the greedy algorithm in Mod- 
elFinder (Kalyaanamoorthy et al. 2017) with the Bayesian 
information criterion and the FreeRate heterogeneity. We 
set the perturbation strength (p) and the number of itera- 
tions since the last best tree was found (c) to 0.3 and 1000, 
respectively. The SH-aLRT (Guindon et al. 2010) and the 
ultrafast bootstrap (UFBoot; Minh et al. 2013) were used 
to assess the support of the branching patterns estimated 
in the phylogeny with 0.99 of the minimum correlation 
coefficient and 1,000 of the maximum number of itera- 
tions. Bayesian analyses were performed with MrBayes 
3.2.1. The best-fitting substitution models were selected 
by jModelTest under the Bayesian information criterion 
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(Suppl. material 1: table S3; Posada 2008). The Markov 
Chain Monte Carlo (MCMC) chain was run for 150 mil- 
lion generations using parameters unlinked among parti- 
tions and sampled every 100 generations. We used Tracer 
v1.5 (Rambaut and Drummond 2009) to monitor the mix- 
ing of the MCMC chains. A burn-in sample of 375,000 
trees was discarded, and a 50% majority rule consensus 
tree was computed with the remaining trees. 


Divergence time estimation 


An uncorrelated lognormal relaxed molecular clock mod- 
el was used to estimate divergence time in BEAST v1.8.1 
(Drummond and Rambaut 2007). The birth-death specia- 
tion process was chosen as the tree prior. Partitioned strat- 
egies (Brandley et al. 2005) were incorporated in BEAST 
analyses, and each gene was used as a separate partition. 
For each partition, the specific model of evolution was 
recommended by jModelTest (GTR+I+G for 12S, 16S, 
NADHI1, and COI, HKY+G for 18S, ITS2, and H3). We 
ran the MCMC for 50 million generations and sampled 
every 1000 generations. The maximum clade credibility 
tree was computed using TreeAnnotator v1.8.0 based on 
the remaining trees, after discarding the first 25% of the 
yielded trees as burn-in. Tracer v1.5 was used to determine 
convergence and measure the effective sample size (>200) 
for all parameters. We dated the tree of all globally sam- 
pled Pardosa species. For the molecular clock analysis, we 
used the minimum ages based on fossils of Lycosidae (15 
Ma; Iturralde-Vinent and MacPhee 1996; Penney 2001), 
Oxyopidae (43 Ma; Wunderlich 2004b; Magalhaes et al. 
2020), Thomisidae (43 Ma; Wunderlich 2004c; Magalhaes 
et al. 2020), and Selenopidae (53 Ma; Penney 2006; 
Magalhaes et al. 2020) as calibration points (Suppl. materi- 
al 1: table S4; Renner 2005; Donoghue and Benton 2007). 
The Agelenidae, Thomisidae, Oxyopidae, Psechridae, Tre- 
chaleidae, Selenopidae, and some Lycosidae species were 
used as outgroups. Their gene sequences were available 
from our study and GenBank (Suppl. material 1: table S1). 


Biogeographical reconstruction 


The biogeographical history of Pardosa was reconstructed 
in RASP v3.0 (Yu et al. 2010, 2015) using the Bayesian 
binary MCMC analysis (BBM; Sanmartin et al. 2001), the 
statistical dispersal-vicariance analysis (S-DIVA; Yu et al. 
2010), and the dispersal-extinction-cladogenesis analysis 
(DEC; Ree and Smith 2008). We set two areas for the max- 
imum number of ancestral states at each node. The fos- 
sil-calibrated trees obtained from BEAST analysis, from 
which the outgroups were trimmed, were used for biogeo- 
graphic reconstruction. The distributional data of Pardosa 
Species were available from previous literature, GenBank, 
and information from our samples. Based on the current 
distribution and the phylogeny of Pardosa as well as geo- 
graphic divisions and climate, we defined five geographic 


zse.pensoft.net 


1290 


areas occupied by the genus: southern Africa, north-east 
Africa, West Asiat+Europe/Palearctic, southern and eastern 
(SE) Asia/Orient, and North America/Nearctic. 


Results 
Phylogeny of Pardosa 


The concatenated data set was composed of 2694 nu- 


cleotides (12S: 237 bp; 16S: 309 b 
611 bp; H3: 383 bp; ITS2: 195 bp; 
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with 133 (125 described and 8 as yet undescribed) termi- 
nals within the Pardosa grassland lineages covering most 
of their distribution ranges. Maximum likelihood analy- 


ses indicate that Pardosa consists of four clades (Fig. 2; 
Suppl. material 1: fig. S1). Clades I and HI only include 
P. crassipalpis Purcell, 1903, from South Africa and 
P. sutherlandi (Gravely, 1924) from SE Asia, respective- 
ly; Clade II is composed of the north-east African group 
(P. injucunda (O. Pickard-Cambridge, 1876) and P. nae- 
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Figure 2. Phylogenetic tree of 133 Pardosa species reconstructed using the maximum likelihood method. The numbers at the nodes 
represent bootstrap support values from the maximum likelihood analyses. 
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diuturnal group (DIG) comprise 7 and 2 Pardosa species 
from North America, respectively; the paludicola group 
(PAG) is comprised of 18 species from Europe and North 
America; the alacris group (ALG) includes 11 species 
from SE Asia, Europe, and North America; and the brevi- 
vulva group (BRG) consists of 95 species from SE Asia, 
West Asia, Europe, and North America. Bayesian analy- 
ses suggest that all sampled African species (P. crassipal- 
pis, P. injucunda, and P. naevia) constitute Clade I with a 
low posterior probability (only 0.62), and the phylogenet- 
ic relationships among most species/groups within Clade 
IV are unclear (Suppl. material 1: fig. S2). 


Divergence time 


The fossil-calibrated phylogeny 1s shown in Fig. 3. The 
recent ancestor of Pardosa probably appeared during the 
early Miocene, approximately 19.40 Ma (95% credibili- 
ty interval, CI: 24.72-14.47 Ma). The divergence of the 
South African P. crassipalpis (Clade I) from the species of 
the clades (I-IV) was predicted to have occurred around 
14.18 Ma (95% CI: 18.19-10.74 Ma). Thus, Pardosa 
wolf spiders originated probably at about 19.40-14.18 
Ma (the Middle Miocene Climatic Optimum, MMCO). 
Within Clade IT, the Pardosa species from north-east Af- 
rica diverged from those from SE Asia around 9.64 Ma. 
The split between the SE Asian P. sutherlandi (Clade 
IIT) and the taxa of Clade IV from Asia, Europe, and the 
Americas started about 12.59 Ma. Clade IV first diver- 
sified during the late Miocene, approximately 8.82 Ma 
(95% CI: 10.77—8.73 Ma), and then split into distinct 
species groups during the Miocene. The North American 
diuturna group (DIG) diverged from the nigriceps group 
(PABG) about 8.17 Ma. Diversification of the widespread 
species group (PABG) started around 7.77 Ma. 


Biogeography 


The BBM analysis (Fig. 4) inferred the two possible an- 
cestral ranges occurring during the middle Miocene for 
Pardosa: (1) southern Africa (the marginal probabili- 
ty, MP: 0.6616); and (2) SE Asia/Orient (MP: 0.2391). 
The inference from S-DIVA (Suppl. material 1: fig. S3a) 
favored southern AfricatSE Asia (MP: 0.6534), SE 
Asiatnorth-east Africa (MP: 0.2066), and southern Af- 
rica+North America/Nearctic (MP: 0.0586) as the most 
likely ancestral areas of the genus, whereas southern Af- 
rica+SE Asia was preferred under DEC (Suppl. material 
1: fig. S3b). The evolutionary routes of Pardosa inferred 
from the BBM (Fig. 4), S-DIVA (Suppl. material 1: 
fig. S3a), and DEC (Suppl. material 1: fig. S3b) were gen- 
erally concordant. The ancestral reconstructions suggest 
that SE Asia is an important dispersal center for Pardosa, 
with dispersal from SE Asia to North America, north- 
east Africa, and Europe at around 10.59, 9.64, and 5.28 
Ma, respectively. The three intercolonial evens between 
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South-east Asia and North America were observed during 
the period of about 10.59-5.31 Ma. The inter-colony dis- 
persal between Europe and North America occurred be- 
tween 5.48 and 2.80 Ma. 


Discussion 


A detailed exploration of the evolutionary history of Par- 
dosa requires a robust phylogenetic framework (Graybeal 
1998; Zwickl and Hillis 2002). Current knowledge of the 
evolutionary relationships among Pardosa species 1s 
based on assessments of morphological similarity or phy- 
logenetic analyses of lycosids, including very few species 
of Pardosa (e.g., Zehethofer and Sturmbauer 1998; Vink 
et al. 2002; Muster and Berendonk 2006; Murphy et al. 
2006; Piacentini and Ramirez 2019; Roslin et al. 2022). 
This study first comprehensively addresses the phyloge- 
ny of global grassland species of Pardosa using multiple 
nuclear and mitochondrial markers. The topology of the 
maximum likelihood tree (Fig. 2) is similar to that of the 
dated tree (Fig. 3) and provides a well-supported hypothe- 
sis for the relationships among the major clades/groups of 
Pardosa within our dataset. The biogeographic histories 
of Pardosa suggest that the evolution of this genus was 
strongly affected by grassland expansions resulting from 
historic climatic and environmental shifts. These findings 
mirror other studies highlighting the vital role of the evo- 
lution of grasslands in lycosid biogeography (Jocqué and 
Alderweireldt 2005; Piacentini and Ramirez 2019). 
Lycosidae appeared about 33.80 Ma, after the Eo- 
cene-Oligocene extinction event, but well before the 
grassland expansions. Pardosa (Pardosinae) diverged 
from the clade (Lycosat+Trochosa) around 19.40 Ma, 
which is consistent with the age (25-16 Ma) obtained 
in the phylogenetic analysis of Piacentini and Ramirez 
(2019). The initial diversification of Pardosa spiders 
was about 14.18 Ma, coinciding with the Middle Mio- 
cene Climatic Optimum (17—14 Ma). The time (19.40— 
14.18 Ma) when Pardosa appeared is much later than 
the original age (55 Ma; Stromberg 2011) of C4 grasses 
with the hot and dry tolerance (Linder et al. 2018), and 
the age (30 Ma) when in colder climates, C3 open-hab- 
itat grasses developed the tolerance needed to survive 
frosts (Schubert et al. 2019), but coincides with the 
time (the late Cenozoic) when open-habitat grasses be- 
gan to become ecologically dominant (Stromberg and 
Staver 2022). Four major clades of Pardosa emerged 
during the middle to late Miocene, between 14.81 and 
8.82 Ma. This pattern coincides with the rapid expan- 
sion of open-habitat grasses and the retreat of tropical 
forests due to the cooling and aridification climate and 
the low CO, conditions, as well as fire and herbivory 
activities after the MMCO (Zachos et al. 2001; Strom- 
berg and Staver 2022). Our results support the idea that 
Pardosa is a young lineage of wolf spiders that became 
vagrant and diversified during the Miocene, when the 
grasslands expanded. In the past 8 Ma, including the 
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Figure 3. A dated, fossil-calibrated phylogeny of Pardosa species estimated in BEAST. Yellow dots mark the four fossil calibration 
nodes. “*” indicates stable branches with Bayesian support >0.90. Ma, millions of years ago. Qu, Quaternary; Pl, Pliocene; Mi, 
Miocene; Ol, Oligocene; Eo, Eocene; Pa, Palaeocene; Cr, Cretaceous; Ju, Jurassic. 


global Ice Age (2.6 Ma onward), climatic cooling fur- 
ther accelerated, and grasslands extended as the forests 
retreated and the open and arid environments expand- 
ed. At this time, Pardosa spiders rapidly diversified, 
currently representing the most abundant lycosid genus 
and being dominant in grasslands and some open envi- 
ronments. 

Our ancestral reconstructions suggested that the Par- 
dosa spiders originated in southern Africa or SE Asia/ 
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Orient during the middle Miocene, likely evolving con- 
currently with the spread of grasses and the diversifica- 
tion of herbivores (~17 Ma; Charles-Dominique et al. 
2016; Stromberg and Staver 2022), followed by the geo- 
graphic expansions from middle altitudes of SE Asia to 
high altitudes, and then to North America via the Bering 
Strait, as well as to north-east Africa southwards via the 
Arabian Peninsula, and to Europe westwards via West 
Asia between about 10.59 and 5.28 Ma. This dispersal 
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Figure 4. Biogeography of Pardosa and potential global dispersion routes (arrows). Reconstruction using Bayesian binary MCMC 
(BBM) in RASP v3.0. The colors of pie wedges at each node represent geographical areas inferred to have been occupied by an- 
cestral taxa. Pink circles in the maps indicate the possible ancestral ranges occurring during the middle Miocene for Pardosa. The 


numbers at the nodes represent support. 


pattern is likely associated with the continent-specific tra- 
jectories of grasses (Karp et al. 2021; Kukla et al. 2022). 
For example, the subtropical C3 open-habitat grasses 
first spread to the colder regions (Kukla et al. 2022). In 
the time (about 10-6 Ma), the tropical open-habitat C4 
grasses expanded to form grasslands and savannas at 
low to mid-latitudes (Karp et al. 2021; Lu et al. 2020), 
and the frost-tolerant grasses spread to higher latitudes. 


During around 10.59 to 5.31 Ma, multiple exchanges (at 
least three times) for the spiders occurred between North 
America and South-east Asia. The Pardosa exchanges 
also appeared between Europe and North America during 
around 5.48 to 2.80 Ma. These results indicated that in- 
tercolony dispersals were associated with grassland ex- 
tensions and retreats resulting from climatic shifts such 


as historical glaciations. 
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The strong dispersal capacity and adaptability to 
complex climates and disturbed environments (Richter 
1970; Greenstone 1982; Samu and Szinetar 2002; Joc- 
qué and Alderweireldt 2005; Woolley et al. 2016) led to 
a world-wide distribution of Pardosa spiders that were 
dominant in open habitats within only about 15 Ma. Fre- 
quent ballooning and cursorial dispersal confer their high 
mobility and area expansion (e.g., Richter 1970; Green- 
stone 1982). Our findings indicate that the long-distance 
overwater dispersal and subsequent range expansions of 
Pardosa spiders occurred between South-east Asia and 
southern Africa and between North America and Europe. 
Furthermore, the marked preference of Pardosa wolf spi- 
ders for disturbed habitats, such as clearings in grasslands 
or forested areas (Jocqué and Alderweireldt 2005), made 
them especially successful when arid and open habitats, 
including grasslands, extended over the globe during the 
late Cenozoic. 


Future directions 


This study generated hypotheses regarding the origin 
and dispersal of Pardosa grassland lineages and sug- 
gested that grassland expansions drove its global di- 
versification during the late Cenozoic using the nuclear 
18S, ITS2, and H3 and mitochondrial 12S, 16S, COZ, 
and NADH] loci. Sampling is the process of choosing 
a subset of a target lineage that will serve as its rep- 
resentative. In this paper, the total sampling specimens 
for grassland species and some regions are on the low 
side. However, we hope our study can aid in strategic 
resampling, reflecting known lineage divergences from 
grasslands. Moreover, the taxonomic revisions of Par- 
dosa need to be made in advance. 
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